Under consideration for publication in J. Fluid Mech. 



1 



Swinging and tumbling of elastic capsules in 

shear flow 

S. KESSLER, R. FINKEN, and U. SEIFERT 

II. Institut fiir Theoretische Physik, Universitat Stuttgart, 70550 Stuttgart, Germany 

(Received 27 June 2007) 

The deformation of an elastic micro-capsule in an infinite shear flow is studied numer- 
ically using a spectral method. The shape of the capsule and the hydrodynamic flow 
field are expanded into smooth basis functions. Analytic expressions for the derivative 
of the basis functions permit the evaluation of elastic and hydrodynamic stresses and 
bending forces at specified grid points in the membrane. Compared to methods employ- 
ing a triangulation scheme, this method has the advantage that the resulting capsule 
shapes are automatically smooth, and few modes are needed to describe the deformation 
accurately. Computations are performed for capsules both with spherical and ellipsoidal 
unstressed reference shape. Results for small deformations of initially spherical capsules 
coincide with analytic predictions. For initially ellipsoidal capsules, recent approximative 
theories predict stable oscillations of the tank-treading inclination angle, and a transition 
to tumbling at low shear rate. Both phenomena have also been observed experimentally. 
Using our numerical approach we could reproduce both the oscillations and the tran- 
sition to tumbling. The full phase diagram for varying shear rate and viscosity ratio is 
explored. While the numerically obtained phase diagram qualitatively agrees with the 
theory, intermittent behaviour could not be observed within our simulation time. Our 
results suggest that initial tumbling motion is only transient in this region of the phase 
diagram. 



1. Introduction 

The dynamics of soft objects such as drops, capsules and cells in flow represents a long- 
standing problem in science and engineering, but has received increasing interest recently, 
in particular due to its relevance to biological, medicinal and microfluidic applications. 
This problem is challenging from a theoretical point of view, because the shape of these 
objects is not given a priori, but determined dynamically from a balance of interfacial 
forces with fluid stresses. Improved experimental methods have revealed intriguing new 
dynamical shape transitions due to the presence of shear flow. The phenomenology of 
the dynamical behaviour depends distinctively on the specific soft object immersed into 
the flow with fluid bilayer vesicles and elastic microcapsules as prominent classes. 

Fluid bilayer vesicles assume a stationary tank-treading shape i n linear shear flow , 
if there is no viscosity contrast between interior and exterior fluid (jKraus et aLl 119961 ). 
If the interior fluid or the membrane becomes more viscous, a transition to a tum- 
bling state can occ ur en fc Misbahll2003l: iBeaucourt et aZ.ll20046l : iRioual et aZ.ll2004t 
Noguchi fc Gompp er 2004. 2005 HVlahovska fc Graciall2007l). Tank-treading was obse rved 
experimentally in infinite shear flow ( de Haas et a^.lll997t Kantsler fc Steinberg||2005 ) and 
for vesicles interactin g with a rigid wall ( Lorz et all ^200rt ^Abkar ian et al. 2002[) . where 
a dynamic lift occurs (|Seifertiri999ft ICantat fc Misbah.1999. : .Sukumaran fc SeiferdlioOll 
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Beaucourt et al. 2004 al ). The tank-treading to tumbling transition was obs erved for the 



first time convincingly in a recent experiment ( Kantsler fc Steinberg 20061 ). In addition 
to the tank-treading to tumtjling t ransition, a vacillating or breat hing motion was pre- 
dicte d theoretically (jMisba h '2006") and observed expe rimentally (jKantsler fc Steinberg! 
2006h and in simulations ( |Noguchi fc Gomppeii 120071 ) The theoretical descripti on has 



been expanded recently beyond first order in the shear rate (jLebedev et aLl 120071 ). 

In contrast to fluid vesicles, microcapsules exhibit a finite shear elasticity, since their 
membrane is che mically or physical ly cross-linked. This includes both artificial poly- 
merised capsules ( Walter et aZ.|[200ll ) and red blood cells (RBCs), whose membrane is 
compose d of an i ncompressible lipid bilayer underlined by a thin elastic cytoskeleton 
|Mohandas fc Evans .1994 ) . The resistance to shear leads to qualitatively different be- 
haviour, such as preventing the prolate to ob late shape transition in viscous fluid vesicles 
in channel flow (jNoguchi fc Gomppen 120051 ). Perhaps most surprisingly, it also leads to 
qualitatively diffe rent instabilities lik e wrinkling first observed experimentally on poly- 
merised capsules (jWalter et aLl 120011 ). which has t o be distinguished fro m the transient 



creasing formation observed later on fluid vesicles (jKantsler et aLll2007T ). The formation 



of the short length scale wrinkles is driven by compressive stress imposed on the mem- 
brane by the flow, while the selectio n of the short length sca le is due to a balance between 
elastic stresses and bending forces ([Finken fc Seifer d l2006h . 



When the unstressed initial shape of the cell is not spherical, material elements of the 
membrane are deformed when d is place d from their initial position. This shape memory, 
suggested for RBCs by iFischer ( 2004I ). leads to a oscillation of the inclination angle 
superimposed on the t ank-treading motion and an intermittent regim e between tank- 
treading and tumbling (jSkotheim fc Secombll2007HAbkarian et aLll2007l ). For a review of 
the tank-treading behaviou r of sof t capsules in shear flow, we refer the reader to the first 
two chapters of IPozriki dis (j2003ah . 

Analytic descriptions of the rather complex motion of capsules and vesicles is only pos- 

sible for asymptotic cases, e.g. in the quasisphcrical limit ('Barthcs-Bicscl'l980':'Barthcs-Bicscl fc RallisonI 
1981; Soifert 1999a; Misbah.2006; Finken fc Seifcrt.2006 : Lobedov et a;..2007; Vlahovska fc Gracial 
2007 ) or by restricting the numb er of degrees of freedom ( Keller fc Skalakll982t " Rioual et al\ 



2004 Iskotheim fc Secombll2007t ). One therefore has to reso't to numerical methods for 



more complex geometries. 

For the dynamics of vesicles existing solvers which treat the flo w at a continuum level 



either em ploy a discrete triangulation scheme (jKraus et al. 



(iBiben fc Misbah.2003i ). An alternative route was taken by Noguchi fc Gompper ( 20041 . 
20051) . where the membrane is dynamically triangulated and the flow is modelled by 



19961) or phase field niodels 



discrete effective fluid particles. 

N umerical simul a tions of capsules were firs t performed in an axisymmetric geome - 
trv (iLi et all Il988t iLevrat-Maurin et all 1 19931 : iLevrat-Maurin fc Barthes-Biesell Il994[ ). 
IPozrikidid (|l995h developed a method for simulating three-dimensional deformations of 
initially spherical capsules i n shear flow using a boundary element formulation. This 
method was later refined bv lRamanuian fc Pozrikidi"s ( 19981 ). who also observed oscilla- 
tions of the inclination angle for ellipsoidal capsules. However, their method was plagued 
by numerical instabilities for high and low deformations due to the degradation of the 
grid. Further improvement of the boundary element method allows the stable simulation 
of tank-treading and tumbling motion of highly flattened capsules only by numerically 
smoothing the surface (jPozrikidis 2003 &) . Part of the numerical difficulties might be due 
to the evaluations of bending moments: Calculating the local mean curvature requires tak- 
ing second derivatives of the shape functions, which become inaccurate using finite differ- 
ence schemes. Newer approaches, such as the spectral boundary algorithms proposed for 
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droplets ( Wang fc Dimitrakopoulosll2006h and solid particles ( Pozrikidi^ 2006h . therefore 
use higher order basis functions on the triangulated surface. While these methods are very 
versatile, the details of these methods are rather complex. Sui table interfacial smoothing 
is sti ll needed to ensure numerical stability of the method ( Wang fc Dimitrakopouloa 
20061) . 



It is therefore the purpose of this paper to augment these approaches with a global 
spectral method. In thi s method th e shape of the capsule is expanded into a set of 
smooth basis functions (jBovdl 120011 ). This has the advantage that the resulting shape 
is automatically globally smooth, which reduces the discretisation error especially in 
higher derivatives, such as the local mean curvature. Since the derivatives of the basis 
functions are analytically known, it is easy to evaluate the elastic tensions and bending 
moments at a grid of collocation points. These marker points are material points of the 
underlying connected membrane. Rather than treating the hydrodynamics in a boundary 
layer formulation, we expand the Stokes flow similarly in terms of smooth basis func- 
tions. Requiring force balance at the collocation points yields the equation of motion 
of the membrane. This scheme is used to systematically explore the dynamic behaviour 
of capsules in shear flow, focusing on init ially non-spherical capsules as considered in 
the analysis by Skotheim fc SecombI ( 20071 ). Although the overall phase behaviour of the 
capsules is captured qualitatively by their model, we could not observe the predicted 
intermittent behaviour. An analysis of the capsule dynamics suggests that the initial 
tumbling motion is only a transient towards a stable tank-treading motion. 

This paper is organized as follows: In section [2l after introducing notions of differential 
geometry and elasticity, we define the problem rigorously. In section [3l we develop the 
spectral algorithm to calculate the dynamics of an clastic capsule. After extensive testing 
for analytically known limit cases in section 31 we apply our method to ellipsoidal capsule 
in section [5] Our findings are summarised in section [6l In the appendix, we recall the 
relevant differential geometry for deformed capsules. 



2. Problem formulation 

We consider the dynamics of a closed capsule that is embedded in an infinite ambient 
flow with viscosity rj° (see figure [ij. The clastic membrane encloses a second fluid with 
a different viscosity r/', defining the dimcnsionlcss viscosity contrast 

e^rf/rf. (2.1) 

In the absence of the capsule we assume a prescribed external flow u°°{x). Because of its 
small thickness we consider the membrane as a two dimensional interface that separates 
the two fluids. On the typical length scales considered inertial effects of the membrane 
are negligible. 

Strain and Curvature 

In order to describe the two dimensional membrane, which is embedde d in three di- 
mension al space, we reca ll some important terms of differential geometry ( Frankell 2004 : 



iNIarsde n fc HugheslflOsl . For mathematical details of the quantities used here, we refer 



the reader to the appendix. A comprehensive summary of i nterfacial propert ies in the 
context of membranes in hydrodynamic flow can be found in IPozrikidisI (|200l[ ). 



Since we consider closed membranes with the topology of a sphere S'^, we can label 
the material points of some reference membrane by spherical coordinates (i?, i^a). Note 
that for an arbitrarily deformed membrane the material point labelled by the Lagrange 
coordinates {•&, tp) will be moved, and {d, tp) will not remain spherical coordinates. 
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The location of the membrane M at thiie t is given by the shape function r{'d,ip;t). 
Length and angles on the membrane are measured by the first fundamental form or 
metric tensor g with covariant components gij (|A 3p . The second fundamental form or 
extrinsic curvature tensor k with covariant components fcy (|A Sp measures how the unit 
normal vector of the surface changes its direction, when one moves along the membrane. 
The mean curvature H is defined as the arithmetic mean of the principle curvatures, 
which are both the eigenvalues of the curvature tensor k and the inverse of the principle 
curvature radii (|A 9[) . In our convention jASj, the mean curvature H = 2/r of a sphere 
with radius r is positive. First and second fundamental forms completely fix the shape of 
the given membrane and therefore contain all information about the membrane shape. 

In order to describe a deformation and an elastic response, we have to specify an 
unstressed reference membrane Ai^ci given by the shape function R{'d,ip). The corre- 
sponding metric is denoted by G and defined analogously to the metric tensor g. The 
Lagrangian strain tensor e is in covariant components given by half the difference of 
metric gij and reference metric Gij (|A 10|) a nd measures the change in length elements of 
the membrane upon deformation (appendix, Marsden fc Hughe^ll983l ) . The strain tensor 
£ holds all information about the deformation and will be used to define an elastic energy 
density. 

Finally, the surface dilation J (jA measures how an infinitesimal patch of area dA 
on the reference membrane is changed upon deformation into the patch of area da (jA 7|) 
and can be expressed by the ratio of determinants of metric and reference metric (jA Sp . 



Constitutive laws, energy, force, stress 

The deformation of the membrane from its unstressed shape costs energy, which can be 
quantified by the underlying constitutive law. In general we consider resistance against 
shear, dilation, and bending. Several elastic models for thi n shells and membra n es are 



considered in the literature. A short overview is found in iBarthes-Biesel et all (|2002h 



and in the first two chapters of Pozrikidii ( 2003q ) . These references directly connect the 
deformation with stresses and bending moments. We prefer to de rive the constitutive law 
from the elastic energy, as outlined in Marsden fc Hughes ( 1983[ ). 

Deformation from the reference shape -R(i?, f) to a shape r{'0, ip; t) costs energy 7i[r]. 
We only consider constitutive laws derived from an energy density h[r\, i.e. 



n[r]^ J dAh[r] 



(2.2) 



Variation of the shape r by 5r while leaving the reference shape fixed induces a variation 
of the total free energy 



Sn= I dA ■5r = - f daf -Sr, 



(2.3) 



which defines the elastic surface force density / on the membrane by a functional deriva- 
tive 



(2.4) 



J 5r 

We now specialise to deformation energies, which can be written as the sum of a purely 
elastic and a bending part, 

n[r]^n^,[r]+n^[r]. (2.5) 
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To illustrate the method, we use the specific constitutive law for the elastic free energy 

7^01 = y dA^ ^"^^^ (tr£)V2/idet£^ , (2.6) 



and the curvature term (jHelfrichlll973f ) 

I da^{2H-Cof ^ I dAJ^{2H-Co)\ (2.7) 

for the bending energy. Here, k is the bending rigidy and Cq is the spontaneous curvature, 
while A and /i are the 2d-Lame coefficients in Hooke's law valid for small deformations. 
These coefficients correspond to the surface shear modulus Gs = ti a n d the surface 
Poisson ratio = A/(A + 2/j,) in the notation of iBarthes-Biesel et all (|2002l ). Other 
constitutive laws can trivially be implemented. 

Hydrodynamics 

For all experimental setups the Reynolds number is small and the flow is governed by 
the linear Stokes equations. The velocity it" (a;) of the inner (a = i) and outer (a = o) 
fluid at the position x is determined by the incompressibility relation 

V-it" = 0, (2.8) 

the linear momentum equation 

- Vp" +77"Am" = (2.9) 

with the isotropic pressure p", and by appropriate boundary conditions far away from 
the capsule and at the membrane. 

The flow has to be regular everywhere and continuous across the membrane when 
assuming a no-slip boundary condition. The jump in hydrodynamic traction between 
both fluids is compensated by the elastic forces at the membrane 

f{d, ^) ^ [T'(r(7?, ^)) - T°(r(^, ^))] • n{d, ^) , (2.10) 

where T" is the inner and outer hydrodynamic stress tensor with Cartesian components 

^ -%p" + r {d.u^ + dX) ■ (2.11) 

Since assuming no-slip boundary conditions, the velocity is continuous across the mem- 
brane, and the membrane is advected with the flow 

Far away from the capsule the outer flow coincides with the undisturbed external flow 

u°{x)^u°°{x) iov\x\^oo. (2.13) 

Since the Stokes equations are linear, we can split the total flow into a sum of the 
undisturbed flow and an induced flow 



(2.14) 



where the homogeneous boundary condition uf^-^^{x) — > far away from the capsule is 
easy to implement. 

For specific applications, we employ a linear shear flow (figure [1]) 

w°°(x)=7j/e, (2.15) 
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with shear rate 7. The equations of motion of the membrane are fuUy determined by 
Stokes' equations ()2.81 12. 9|) . the force balance with the elastic forces ()2.10|) . and the 
boundary conditions (j2.131 12.12p which include the membrane advection (right hand side 

Dimensionless parameters 

The motion of the capsule is governed by a number of dimensionless parameters: The 
volume V of the capsule remains constant and defines a length scale Rq 

V^fRl, (2.16) 

which will be used as the unit length. The elastic energy density is given by the elastic 
moduli depending on the given constitutive law. In our case we use the shear elasticity 
fj, to define an energy scale ^Rq . The remaining elastic parameters can thus be cast in a 
non-dimensional form by defining the two dimensional Poisson number 

A 

1/ = . 

X + 2n 

and the non-dimensional bending rigidy 



(2.17) 



and spontaneous curvature 



(2.18) 



Co = ^Co. (2.19) 
The viscosity r]° can be used to define a time scale R^if / fi, giving the capillary number 

X^^i- (2.20) 
Finally, the viscosity contrast e has already been defined in (|2.ip . 



3. Spectral Method 

We now develop a method to numerically solve the nonlinear equations of motion. 

Spectral method 

To transform the shape function to spectral space, we expand its Cartesian compo nents 
r^('d, i p) = r{d, (fi) ■ e,; into sca lar spherical harmonics Yi^{'d, (f) with Z ^ 0, |m| < Z (jRose 
1957t Brink fc Satchler 1968 ). More generally, we consider the spectral expansion of a 



scalar function / 

/(^,^) = ^/'™V'(^,^). (3.1) 

Im 

Since the set of all spherical harmonics form a complete and orthonormal basis on the 
sphere S"^: 



duj y;"(i9, (/p)l^r i^, ip) = S„,„,,Su' , (3.2) 
the spectral coefficients are in principle obtained by the integral 

/ dojYnd,^)f{d,^). (3.3) 
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Here dui = sin'dd'ddip denotes the area element on the sphere S*^ and the superscript 
star indicates the complex conjugate. For numerical applications, however, this integral 
must be replaced by a discrete sum. For the hydrodynamic part it will be quite useful to 
extend the spherical harmonics to / < by defining y^^^j^^ = Y^K 

A function whose spect ral coefhcients fim vanish for I ^ b is called a bandlimited 
function with bandlimit b (IHealv et a/.l 120031) . We solve the dynamics of the membrane 
by restricting to the space of bandlimited shap e functions. Since the spectral amplitude 
of smooth functions decays exponentially with I (|Bovdll200lh . this scheme is very accurate 
already for low b. 

To get the expansion coefficients out of a given bandlimited function, we choose a 
finite number of collocation markers (iJj, i ~ at the membrane. A scalar 

function ip) living on the membrane is then determined by its values at these points 
ji = f{'3i,ipi), whereas in spectral space this function is represented by the spectral 
coefficients /'™ up to the bandwidth b. Transformation from spectral space to real space 
is easily done by evaluating the spherical harmonics at the collocation points 



r \ ^ rlrn\^m / n \ \ ^ Alrn rim 



(3.4) 



hn 



hn 



If the transformation matrix is square and regular, the inverse transformation can 
be obtained by 



(3.5) 



where A'"^ is t he inverse of A. H owever, often more collocation points than spectral 
modes are used ([Healv et aZ.ll2003l ): A natural choice is to define the collocation markers 
i = (j, k) equidistantly in "d and i^, where we shift t he values "d a way from the poles to 
avoid numerical problems in the vicinity of the pole (|Bovdll200l[) 



V(j,fc) 



_ (2j + l)7r 



4& 



kn 



(3.6) 
(3.7) 



with j, = 0, . . . , 26 — 1. With this choice the number of markers n = 46^ is larger than 
the number of spectral modes b^. In this case has to be replaced by the Moore- 
Penrosc-Pscudoinverse ( Swarztrauber fc Spotj 2004 ). 

The expressions of the metric and curvature involve first order and second order deriva- 
tives, respectively. The main advantage of spectral methods is that t he de rivatives of the 
basis functions arc known algebraically (Rose 1957; Brink & Satclilcrl ll968 ). and therefore 
differentiation can be performed with high accuracy for bandlimited functions. Similarly 
the integral over the function / is evaluated numerically via 



(3.8) 



dLuf{d,ip)^Vi7Tfl 



which follows easily from Yi^{-d, ip) = I/V^tt. Once the derivatives of the shape functions 
are known, all further geometrical computations are performed at the collocation points 
in physical space. It is straightforward to numerically calculate the energy density of a 
given shape for a given constitutive law with high accuracy. Similarly, the variation of 
the energy density is evaluated for a given shape and variation of the shape function. 
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Elastic forces 

To get the force density at the collocation points, we need to calculate the variation of 
the total free energy for special variations of the shape. The spectral coefficients of the 
force density are obtained in the most direct way when we choose the specific variations 

sin I? 

6rl^{d,^) = =_y;™(z9,^)e,, I ^ 0, . . . ,b; m ^ . . . ,1 ; i ^ x,y, z . (3.9) 

The variation of the metric and curvature tensor resulting from this shape variation can 
be evaluated easily using the known derivatives of the shape function. Since the deriva- 
tives of the energy density with respect to metric and curvature are known analytically, 
the variation of the elastic energy ()2.3p can be calculated easily.For the specific choice of 
Sfhn P-9P - this yields directly the Cartesian spectral force components 

6-H = - [ daf- 5rl, = f d^ f ■ e,l^""(^, ^) = /L • (3.10) 



Hydrodynamics 

We follow a similar strategy for the hydrodynamic part of the problem. To solve the hy- 
drodynamic equations, we choose a complete set of basis functions in three dime nsional 
space that automatically fulfi ll Stokes' equations. These so called Lamb modes ( LambI 
ll932HHaDDel fc BrenneI^ll983^ can be defined as appropriate linear combinations of vec- 
tor spherical harmonics. Spherical coordinates r, and 4> in the laboratory frame as well 
as the corresponding basis vectors Br , eg and are best suited for calculations concern- 
ing the Lamb modes. We stress again the difference between the Lagrangian coordinates 
{'d,<f), which serve as markers for the material points, and the angles {9,(j)), which are 
spherical coordinates in the laboratory frame. With help of the surface gradient on the 
sphere 



eedg 



smf 



the vector spherical harmonics arc given by (|Morse fc Feshbachlll953f) 



(3.11) 



1 



yi(ITT) 
1 

V^OTT) 



Br X V^l^™( 



(3.12) 

(3.13) 
(3.14) 



The Lamb modes then read 
, 1 



^7L,(^^,0)^-y^(7Tl)*^(^,'/')• 



{i + i)^Ji{i + \)<i 



I V^i' 



,(3 


15) 


(3 


16) 


(3 


17) 



With the scalar spherical harmonics and the Lamb modes as basis functions the pressure 
and velocity fields can be expanded for the external and internal as well as the induced 
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flow 

p"(r,0,^) =77"^P;V'^r(^,^), (3.18) 

Im 

Irn 

Due to the regularity of the induced flow at the origin and the boundary conditions at 
infinity (|2.13p . the sums are restricted to Z for u\^^ and Z < — 2 for 

The remaining boundary conditions (|2.12p and force balance condition (|2.10p result 
in linear equations for the coefficients of the Lamb modes. For our choice of collocation 
points, we have an overdetermined system that can be solved in a least square sense. 
This way we get the hydrodynamic flow for any given deformation of the capsule. 

Since the capsule is simply advected with the flow (|2.12p . we perform Euler steps 
with a sufficently small time step to determine the dynamics of the capsule. Higher 
order methods such as a second order Runge-Kutta have been tested, but did not yield 
significant improvements in terms of stability and overall simulation time. 



4. Spherical Capsules 

To test our method, we compare with quasisphcrical results that can be obtained 
analytically. These tests comprise the relaxation dynamics of a capsule to its spherical 
reference shape in quiescent fluid and the stationary deformation of a capsule in linear 
shear flow for low shear rates. 

Sma ll deformations of an initially spherical capsule relax exponentially in time. lRochal et 
(l2005l) lave identified the normal modes, and calculated the corresponding relaxation 



times. The relaxation modes are linear combinations of vector spherical harmonics and 
are likewise labeled by l,m. For each set of angular momentum numbers, three relax- 
ation normal modes exist, which have been termed "stretching", "bending", and "shear" 
mode, respectiv ely. The correspondin g relaxation times arc obtained from the eigenvalue 



equation (24) of lRochal et al\ (120051 ). 

As a numerical test of our code, a spherical capsule was deformed in the direction of 
a normal mode, and the time constant of the subsequent relaxation to equilibrium was 
extracted. In figure [H we compare the relaxation times of selected modes as a function 
of the Lame parameter A with the theoretical predictions, with excellent agreement. 

Switching on the shear flow with a low shear rate, the capsule relaxes into a stationary 
shape with a tank-treading motion. The station ary deformations have been calc ulated 



to first and second order in the deformation bv iBarthes-Biesel fc RallisonI (|l98lh . The 
deformation of the capsule is measured by the time-dependent Taylor deformation pa- 
rameter 

where L is the longest and S is the shortest distance of the membrane from the center 
(see figured]). In the long time limit, D assumes a stationary value Dq, which is to first 
order in x (Co = 1) 



(4m + 3X)RoV° 
4i/+1 + 2k(i^ + 5)^ + X)fi + j^{3X + 5ij) 



^"~7 .. , 1 , or./., , , , 2k /o n . (4-2) 



In our simulations, the unstressed capsule is subjected to an abruptly started linear 
shear flow. Numerically, the deformation and inclination are not calculated directly ac- 



10 S. Kessler, R. Finken and U. Seifert 

cording to definition (j4.ip . As described bv iRamanujan fc Pozrikidi^ ( 1998[ ). we find it 
more suitable to use instead the deformation parameter of an ellipsoid with the same 
tensor of inertia. In the small deformation limit, both definitions arc equivalent. 

The simulations were performed by using a bandwidth of 6 = 11 corresponding to 
121 modes, and a total of = 46^ — 484 grid points. The bending energy was chosen 
to prevent the formation of wrinkling modes within the bandwidth. The time step was 
chosen small enough to yield a stable evolution. Typical values were At ^ 1/(10007). The 
time evolution of the deformation parameter of an initially spherical capsule is shown in 
figure [3] for different shear rates and fixed elastic parameters. To illustrate the stationary 
tank-treading motion of the capsule, we also show the distance of a marker point on the 
membrane to the center of the capsule as a function of time. This distance oscillates 
with twice the tank-treading frequency between L and S. At low shear rate, we observe 
a monotonic relaxation of the deformation to its stationary value Dp, while for large 
X and e a pronounced over-shoot is observed (c.f. Ramanuian fc Pozrikidii 1998 ). The 
numerical deformation Dq clearly follows the asymptotic prediction (|4.2p for low shear 
rate. For large shear rates, the simulation data show deviations from linear behaviour. 
These are likely due to the rotational part of the linear shear flow, since for non-zero 
vorticity the stationary deformation is a non-linear function of the shear rate even in the 
quasispherical limit. 



5. Ellipsoidal Capsules 

After successfully testing our spectral method by means of a spherical initial or refer- 
ence shape, we can continue to investigate the dynamics of capsules with an ellipsoidal 
initial shape. This case is both experimentally more realistic, since synthcsised capsules 
arc never perfectly spherical, and leads to a richer dynamical behaviour of the capsules. 
In a non-spherical reference shape, the membrane points arc no longer equivalent to each 
other. During the course of a tank-treading motion, a membrane clement is therefore 
periodically deformed, which costs elastic energy. Any membrane clement therefore en- 
ergetically prefers its initial position (or one of the equivalent positions by symmetry of 
the reference shape). This effect is calle d "s hape me mory" , and plays an im portant role 
also in the dynamics of red blood cells |Fisch cr 2004; Watanabe et aL|[2006[) . 

The shape memory effect is the cause of oscillations of the deform ation and the incli- 
nation angle in the tank-tre ading state, as obs erved experimentally bv lChang fc Olbricht 
(Il993l). lwalter et all (l200ll). and lAbkarian et al. (200n ). and also found in the simula- 
tions by Ramanuian fc Pozrikidid ( 1998 ). Recentlv. a modified Keller-Skalak type theory 
was proposed ( Skotheim fc Secomb 20071 ). which explains this behaviour qualitatively. 
Their model also predicts an oscillating tank-trcading motion at large shear rate, and 
a tumbling motion at lower shear rate. In the tumbling regime, a tracer particle on the 
membrane oscillates around a fixed position with respect to the capsule shape. In the 
intermediate shear rate regime, intermittent motion, which alternates between tumbling 
and tank-treading, is predicted. Although direc t experimental evidenc e for this behaviour 
is missing, indirect evidence was provided by lAbka rian et ali (|2007l ). who discovered a 
hysteresis of the transition from tumbling to tank-treading and the reverse transition by 
increasing or decreasing the shear rate, respectively. We use our spectral method to sys- 
tematically explore the full phase diagram in a large range in shear ra te as well as viscosity 
contrast. Thus the quantitati ve accuracy of the reduced mod el in ISkotheim fc Secomb 
( 2007 ) can be tested. While Ramanuian fc Pozrikidid (1998) observed the onset of a 
tumbling motion for low sphericity, due to the formation of cusp-like instabilities in the 
shape the simulations never went beyond half a tumbling motion. Grid distortion also 
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requir ed the use of exphcit numerical smoothing in more recent simulations ( Pozrikidi^ 



200361 ). Since bending rigidity is included in our method, the formation of cusps can 
be suppressed, leading to a more stable algorithm. The cutoff at a finite bandwidth in 
our algorithm also effectively implements numerical smoothing. Nevertheless high order 
modes accumulate numerical errors during the simulation run in particular at large shear 
rates, thereby limiting the maximum simulation time. 

Phase Diagram 

Our numerical results for the overall phase diagram are summarised in figure^ where the 
dynamical behaviour is plotted as a function of the inverse dimensionless shear rate 
and the viscosity contrast e. At low shear rate, the hydrodynamic forces are too small 
to overcome the energy barrier present for a tank-treading motion due to the shape 
memory effect. Therefore, capsules tumble at low while an oscillating tank-treading 
behaviour is stable at large x- We also observe a transient dynamics from tumbling to 
tank-treading for large viscosity contrast e, which will be discussed below. Although this 
transient dynamics might be taken as indications of intermittent motion, we could not 
find conclusive evidence during the time of our simulation runs. In particular, we never 
observed a transition from tank-treading to tumbling. Also shown in this figure is the 
phase diagram for the analytic model bv ISkotheim fc Secomb ( 2007 '). The qualitative 



agreement, apart from the apparent lack of intermittent behaviour, seems to be rather 
good, given the crude dynamics implemented in the reduced analytical model. Only at 
large viscosity contrast pronounced differences in the shape of the phase diagram start to 
feature. Closer inspection of the data reveals significant oscillations of the axis lengths, 
which are fixed in the reduced model. These bre athing modes may alt e r the intermittent 
character of the capsule motion in the model bv ISkotheim fc SecombI (|2007h . 



For all simulations the deformation remained well within the range of validity of 
Hooke's law. The extension ratios (c.f. appendix [X| never deviated more than 5% from 
unity. The results are therefore not susceptible to the specific elastic law in this regime. 

Oscillation Amplitudes 

We proceed to quantify the oscillations in the tank-treading and tumbling state and 
investigate the transient dynamics below. For the definition of inclination /3 and tank- 
treading angle a see figure [1] They are defined as the angles between the flow direction 
and the maximal radius or a marker point, respectively. Initially, these angles are chosen 
to lie in the interval [0, tt[. To make both a{t) and P{t) continuous functions of time, 
we manually add 2tt after a full rotation, thereby keeping the variation of each angle 
between subsequent time steps as small as possible. As a consequence, the angles can 
assume values outside the interval [0, 2tt[ at later times. This allows to count the number 
of full rotations, which would not be possible if all angles were restricted to [0, 2tt[. The 
phase shift of a material point away from its elastic minimum 

S{t) = {a{t) - m) - HO) - m) (5.1) 

is called phase angle. 

In the tank-treading regime the inclination angle f3 oscillates around a stationary value 
Po with amplitude A/3 while the tank-treading angle a or phase angle S changes monoto- 
nously with time (see figure [5]). In the tumbling regime, the inclination angle (3 changes 
monotonously with time, while the phase angle 5 oscillates around a fixed value Sq with 
an oscillation amplitude AS. Figure [6] shows parametric plots of inclination vs. phase 
angle for a tumbling and a tank-treading motion, as well as a transition from tumbling 
to tank-treading. The arrows indicates the direction of time in this plot. In Figure [6] (c) 
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and (d) , one can see the transition from a tumbling motion to an oscillating tank-treading 
motion near the transition. Despite intensive search, we have not observed the reverse 
behaviour: A tank-treading capsule never started to tumble. This might indicate that 
the initial tumbling motion is only transient. 

For a fixed viscosity contrast e = 10 figure [7] shows the shear rate dependence of the 
mean inclination angle po in the tank-treading regime, of the mean phase angle Sq in the 
tumbling regime, and of the oscillation amplitudes in both the tumbling and the tank- 
treading regime. For low shear rates, in the tumbling regime, the mean phase angle Sq and 
the oscillation amplitude AS of the phase angle are plotted as a function of the shear rate. 
Whereas the mean phase angle depends only weakly on the shear rate, the oscillation 
amplitude of the phase angle increases for increasing shear rate. For low shear rates, this 
amplitude grows approximately linearly with the shear rate. When the amplitude reaches 
approximately 7r/2, the capsule starts to tank-tread. 

For higher shear rates, in the tank-treading regime, the mean inclination angle (3o 
and the oscillation amplitude A/3 of the inclination angle are plotted as a function of x- 
With decreasing shear rate, the oscillation amplitude of the inclination angle increases 
until it reaches approximately 7r/2, where the transition to tumbling takes place. The 
mean inclination angle also increases, but does not reach tt/A. Perhaps surprisingly the 
oscillation amplitude can be larger than the mean inclination angle /3o in the tank- 
treading regime, implying that the inclination angle is negative for short periods of time. 



6. Summary 

During the last few years, the dynamics of elastic capsules in linear shear flow has 
received increasing attention. Theoretical descriptions restricting the capsule deformation 
to a few degrees of freedom predicted a rich dynamical phase diagram, comprising of tank- 
treading, tumbling and an intermittent motion. While recent experiments have found 
a hysteresis in the tank-treading to tumbling transition for varying shear rate, direct 
observations of intermittent behaviour is lacking so far. 

To investigate elastic capsule systems, while maintaining complete control over the 
underlying equations of motion, we implemented a spectral method to numerically solve 
the equations of motion for a capsule. The capsule deformation is expanded into smooth 
basis functions, leading to accurate estimates on the membrane forces. The code is flexible 
and stable enough to permit simulations for a large range of shear rates and viscosity 
contrasts between inner and outer fluids. 

Using this code, we could quantitatively recover the asymptotically known stationary 
deformations of initially spherical capsules for low shear rate. 

Finally we app hed the numerical me t hod t o the ellipsoidal capsule system similar to the 
one discussed bv lSkotheim fc SecombI ( 2007 ). We observe a stable tank-treading motion 
for large shear rate, in which the inclination angle oscillates with twice the tank-trcading 
frequency. At lower shear rate, or higher viscosity contrast, the capsule starts to tumble. 
We systematically explored the capsule dynamics over a wide range of viscosity ratios 
and shear rates. Th e resulting phase diagram i s qualitatively similar to the theoretical 
predictions made bv lSkotheim fc SecombI (|2007h . with the exception of intermittent dy- 
namics: While dynamic transitions from tumbling to a stable tank-treading motion were 
observed, the reverse transition could not be observed. An analysis of the results suggests 
that the tumbling motion is only transient in this part of the phase diagram. 

Much longer simulation runs and a more detailed analysis near the presumed inter- 
mittency region are needed to decide whether intermittent motion is merely an artifact 
of the reduced dynamics employed bv ISkotheim fc SecombI ( 2007 ). Differences from the 
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phase diagram at large viscosity ratios are likely due to the restriction to a fixed capsule 
shape in the reduced model. 

In conclusion, the spectral method developed in this work is a stable and accurate 
complement to existing numerical methods. It allows the systematic exploration of cap- 
sule dynamics over a wide range of material constants. Theoretical predictions of the 
phase diagram of ellipsoidal capsules are qualitatively confirmed, although quantitative 
differences exist, especially for large viscosity ratios. The nature of predicted intermittent 
behaviour warrants further investigation. 

Financial support of the DFG with in the priority programme SPP 1164 "Nano- and 
Microfluidics" is gratefully acknowledged. 



Appendix. Differential geometry for deformed capsules 

We first recall some definitions of d ifferential geometry that can be found in iFrankel 
( 2004 ) and iMarsden fc Hug hei (Il983h . We cover the two dimensional surface with a 
coordinate net as outlined in section [21 where the coordinates {-d, ip) label the material 
points. The location of the surface at time t is given by the shape function r(i?, ip; t). The 
basis vectors, defined by 

Bi = dir, i = 'd,ip, (A 1) 

span the tangent planes and induce the outward pointing normal unit vector 

n = - ^ -e^XBu,. (A 2) 

\ei} X 

They also define the metric tensor g with covariant components 

9ij = ' ) (-^ 3) 

where we use the ordinary three dimensional Euclidian scalar product. The inverse metric 
tensor g^^ with contravariant components is given by 

g''9jk = S'k, (A 4) 

where Einstein's sum convention is employed. The determinant of the metric tensor 

g = detg (A 5) 

is connected to the basis vectors via 

y/g =\ei) X e^l . (A 6) 

The area of an infinitesimal patch with width di? and length dip in Lagrange coordinates 
is given by 

da = yfgdMip . (A 7) 

The curvature tensor k, defined by 

kij = Bi ■ djTi = ~n • diSj , (A 8) 

measures how the normal vector n changes its direction, when one moves along the 
membrane. Its eigenvalues fci_2 are called principle curvatures and are the inverse radii 
of the principle curvature circles. Trace and determinant of the curvature tensor define 
mean H and Gaussian curvature K respectively. They serve as scalar invariants of the 
curvature tensor and can be expressed by the sum and the product of the principle 



14 



S. Kessler, R. Finken and U. Seifert 



curvatures 

2H = trk = g'^kji ^ ki + k2 , K = det k = kik2 ■ (A 9) 

In order to describe a deformation and an elastic response, an unstressed reference 
shape R{'d, if) has to be specified as was mentioned in section[2] The corresponding basis 
vectors E, normal vector N and metric G arc defined analogously. 

The Lagrangian strain tensor s with covariant components 

£^] = l{9v-G^J) , t = d,^ (A 10) 



meas ures the change in length elements of the membrane upon deformation ( Marsden fc Hughes! 
19831) . 



At each point of the reference membrane there are two orthogonal directions called 
principle extension directions for which the deformation is maximal and minimal. The 
corresponding deformed line elements along these directions remain orthogonal and have 
a stretched length given by the so called extension ratios Ai measured in units of the 
undeformed line elements. 

An infinitesimal material patch on the reference shape with area dA = VGdMip will 
deform into the area element da — y/gd'dd(f on the current shape. The surface dilation J 
is therefore given by the product of the extension ratios 

..^^yi^A.A. ,AU) 

All scalar deformation quantities can be expressed by the extension ratios or equivalently 
by the eigenvalues of the Lagrangian strain tensor e^, which measure the difference of the 
extension ratios from unity 

s, ^U\^-l) . (A 12) 
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Figure 1. Elastic capsule in hydrodynamic flow. — The viscosity of the outer flow and the inner 
fluid are rj" and 77' respectively. Long and short axes of the deformed capsule are denoted by L 
and S. The inclination angle /3 measures the angle between the direction of maximal elongation 
and that of the shear flow = 'yyex- The angle defined by a tracer particle on the membrane 
compared to the flow direction is denoted by a. 
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Figure 2. Comparison of numerically and analytically obtained relaxation time. — The plot 
shows numerically obtained values of the scaled inverse relaxation time rfRa/^r of bending 
modes as a function of the ratio A/^ for different harmonic indexes 1 = 2, Z = 3, Z = 4 an d k = 0. 
The curves are the analytic solution of the secular equation (24) of lRochal et al\ (|2005l '). 
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Figure 3. Dynamics of an initially spherical capsule in shear flow. — (a) and (b) Time evolution 
of maximal, minimal radius and radius of a tracer particle for abruptly starting shear flow at 
time t — for different shear rates x ~ 0.01 (a), x = 0.3 (b) — (c) Time evolution of deformation 
parameter D for different shear rates: x = 0.03 (continuous), x = 0-3 (dotted), x = 3. (dashed 
line). — (d) Plot of the stationary deformation parameter Dq as a function of the dimensionless 
shear rate x- For low shear rates the results of the linear theory f eg. 14.21 full line) are approached 
whereas for higher shear rates deviations are clearly visible. — Constant parameters for (a)-(d): 
viscosity contrast e = 10, Poisson number i/ = 0.5, non-dimensional bending rigidy k = 0.01 
and spontaneous curvature Co = 1. 
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Figure 4. Phase diagram of an elastic capsule in shear flow with the tumbling and tank-treading 
regimes as a function of the viscosity contrast e and the inverse dimensionless shear rate x~^- — 
The full line is a guide to the eye separating the tank-treading (circles) , tumbling (crosses) and 
transient region (diamonds ) for our simulation. Dashed lines indicate the phase diagram due to 
ISkotheim fc SecombI (|2007D for the same parameter set. In the region between the dashed lines 
intermittent motion is predicted. We have not found conclusive evidence for this kind of motion, 
but rather found transient dynamics from tumbling to tank-treading. The numbers correspond to 
following figures, where parts of the phase diagram are examined closer. Geometrical parameters: 
a2 = 03 = 0.9ai, elastic parameters; u = 0.333, k = 0.01, Co = 1. 
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Figure 5. Definition of oscillation amplitudes of phase and inclination angle in tumbling and 
tank-treading states of motion. — (a) shows the tumbling state: Here the inclination angle /3 
changes monotonously while the phase angle <5 oscillates with amplitude A(5 as is indicated by a 
tracer particle. — (b) shows the tank-treading motion, where the inclination angle /3 oscillates, 
while the phase angle 5 and the tank-treading angle q change monotonously. 
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Figure 6. Typical plots of inclination angle (3 vs. phase angle <5 for tank-treading, tumbling 
and transient dynamics. — The site of these plots in the phase diagram are labelled in figure 
[4] by corresponding numbers. — (a) Typical tank-treading motion: /3 oscillates around a stable 
value, while 5 changes monotonously. The arrow denotes the direction in time, viscosity contrast 
e = 13.3, non-dimensional shear rate x = 0.08. — (b) Typical tumbling motion: S oscillates 
around a stable value, while /? changes monotonously, e — 13.3, x ~ 0.025. — (c), (d) Typical 
motions for tumbling to tank-treading transition, e = 23, x ~ 0.2 in (c) and e = 27.5, x ~ 0.2 
in (d) . — The remaining parameters are equal to those used in figure [l] 
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Figure 7. Mean phase angle So (small crosses), mean inclination angle /3o (small filled circles) 
and oscillation amplitudes of phase angle AS (plus) and inclination angle A/3 (circles) for dif- 
ferent shear rates x and a constant viscosity contrast e = 10. — This cut through the phase 
diagram with e = 10 is indicated by the dotted line in figure [l] At low shear rates the capsule 
tumbles with A(5 < 7r/2 where the dashed line indicates a linear behaviour for small x, at higher 
shear rates the capsule tank-treads with A/3 < 7r/2. — remaining parameters as in figure |4l 
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